/*
* INSANE - Interactive Structural Analysis Environment
*
* Copyright (C) 2003-2005
* Universidade Federal de Minas Gerais
* Escola de Engenharia
* Departamento de Engenharia de Estruturas
* 
* Author's email :    insane@dees.ufmg.br
* Author's website :  http://www.dees.ufmg.br/insane
* 
* This program is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License
* as published by the Free Software Foundation; either version 2
* of the License, or any later version.
* 
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
* GNU General Public License for more details.
* 
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
*/
package br.ufmg.dees.insane.materialMedia.constitutiveModel;

import java.io.Serializable;
import java.util.HashMap;

import br.ufmg.dees.insane.analysisModel.AnalysisModel;
import br.ufmg.dees.insane.analysisModel.PolarPlaneStrain;
import br.ufmg.dees.insane.analysisModel.PolarPlaneStress;
import br.ufmg.dees.insane.continuousPointModel.ContinuousPointModel;
import br.ufmg.dees.insane.materialMedia.material.Material;
import br.ufmg.dees.insane.materialMedia.material.PolarVonMises;
import br.ufmg.dees.insane.util.IMatrix;
import br.ufmg.dees.insane.util.IVector;

/** A class representing a PolarVonMisesConstMod
 * @author jamile
 * 16/11/2005
 */
public class PolarVonMisesConstModel implements ConstitutiveModel ,Serializable {
    
	private static final long serialVersionUID = 1L;
	
	/** The tolerance used to verify plastification.*/ 
    private double PVM_TOL = 1.0e-6; 
    
    /** The number maximum of iteractions.*/ 
    @SuppressWarnings("unused")
	private int EP_MAX_ITE = 50; 
    
    /** The plastic multiplier used in associative plasticity.*/
    public final String PLASTIC_MULTIPLIER = "PLASTIC_MULTIPLIER";    
     
    /** The plastic multiplier increment used in associative plasticity.*/
    public final String INC_PLASTIC_MULTIPLIER = "INC_PLASTIC_MULTIPLIER";
    
    /** The derivative of yield surface w.r.t. stress used in associative plasticity.*/
    public final String dFdSIG = "dFdSIG";
    
    /** The plastic deformation used in associative plasticity.*/
    public final String PLASTIC_DEFORM = "PLASTIC_DEFORM";
    
    /**Constructs this Constitutive Model.*/
    public PolarVonMisesConstModel()
    {
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#getLabel()
     */
    public String getLabel() {
        return "PolarVonMisesConstModel";
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#mountCt(br.ufmg.dees.insane.analysisModel.AnalysisModel, br.ufmg.dees.insane.materialMedia.material.Material, java.util.HashMap)
     */
    public IMatrix mountCt(AnalysisModel am, Material mat, HashMap hm) {
    	this.check(am,mat);
        if (Math.abs(((Double)hm.get(this.INC_PLASTIC_MULTIPLIER)).doubleValue()) <= this.PVM_TOL){
            return (am.getDualInternalVariableOperator(mat.getPt()));
        }
        IMatrix Ct = new IMatrix(am.getNumberOfInternalVariables(), am.getNumberOfInternalVariables());
        IMatrix c = am.getDualInternalVariableOperator(mat.getPt());
        IMatrix c1 = c.inverse();
        double coeff = 1.5 * ((Double)hm.get(this.INC_PLASTIC_MULTIPLIER)).doubleValue()/ (((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).YIELDING).doubleValue() + ((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).HARDENING_MODULUS).doubleValue()* ((Double)hm.get(this.PLASTIC_MULTIPLIER)).doubleValue());
        IMatrix P = am.getPMatrix();
        IMatrix H1 = new IMatrix(am.getNumberOfInternalVariables(), am.getNumberOfInternalVariables());
        P.scale(coeff);
        H1.add(c1,P);
        IMatrix H = H1.inverse();
        IVector HdFdSig = H.mul(((IVector)hm.get(this.dFdSIG)));
        double denom = ((IVector)hm.get(this.dFdSIG)).dot(HdFdSig);
        denom += 1.0 / (1.0/((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).HARDENING_MODULUS).doubleValue() - ((Double)hm.get(this.INC_PLASTIC_MULTIPLIER)).doubleValue()/(((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).YIELDING).doubleValue() + ((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).HARDENING_MODULUS).doubleValue()* ((Double)hm.get(this.PLASTIC_MULTIPLIER)).doubleValue()));
        for(int i=0; i<am.getNumberOfInternalVariables(); i++)
            for(int j=0; j<am.getNumberOfInternalVariables(); j++)
                Ct.setElement(i,j,H.getElement(i,j)- (HdFdSig.getElement(i)*HdFdSig.getElement(j)/denom));         
        return Ct;
     }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#mountCs(br.ufmg.dees.insane.analysisModel.AnalysisModel, br.ufmg.dees.insane.materialMedia.material.Material, java.util.HashMap)
     */
    public IMatrix mountCs(AnalysisModel am, Material mat, HashMap hm) {
    	this.check(am,mat);
        return (am.getDualInternalVariableOperator(mat.getPs()));
    }

    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#mountDualInternalVariableVector(br.ufmg.dees.insane.util.IVector, br.ufmg.dees.insane.analysisModel.AnalysisModel, br.ufmg.dees.insane.materialMedia.material.Material, java.util.HashMap)
     */
    public IVector mountDualInternalVariableVector(IVector et, AnalysisModel am, Material mat, HashMap<Object, Object> mp) {
    	this.check(am,mat);    	
    	IVector ep = ((IVector)mp.get(this.PLASTIC_DEFORM));       	
        IVector stress = new IVector(am.getNumberOfDualInternalVariables());
        IVector ee = new IVector(am.getNumberOfDualInternalVariables());
        ee.sub(et,ep);
        IMatrix c = am.getDualInternalVariableOperator(mat.getPs());
        stress.mul(c,ee);
        double f = this.yieldFunction(stress, am, mat, mp);
        if(f > this.PVM_TOL){
            this.yieldGrad(stress, am, mp);
            double denom = ((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).HARDENING_MODULUS).doubleValue();
            for(int i=0; i<am.getNumberOfInternalVariables(); i++)
                for(int j=0; j<am.getNumberOfInternalVariables(); j++)
                    denom += ((IVector)mp.get(this.dFdSIG)).getElement(i)* c.getElement(i,j)* ((IVector)mp.get(this.dFdSIG)).getElement(j);
            mp.put(this.INC_PLASTIC_MULTIPLIER, f / denom);
            mp.put(this.PLASTIC_MULTIPLIER, ((Double)mp.get(this.PLASTIC_MULTIPLIER)).doubleValue() + ((Double)mp.get(this.INC_PLASTIC_MULTIPLIER)).doubleValue());          
            for(int i=0; i<am.getNumberOfInternalVariables(); i++){
                for(int j=0; j<am.getNumberOfInternalVariables(); j++)
                    stress.setElement(i,stress.getElement(i) - ((Double)mp.get(this.INC_PLASTIC_MULTIPLIER)).doubleValue() * c.getElement(i,j)* ((IVector)mp.get(this.dFdSIG)).getElement(j));
                ep.setElement(i,ep.getElement(i)+((Double)mp.get(this.INC_PLASTIC_MULTIPLIER)).doubleValue()*((IVector)mp.get(this.dFdSIG)).getElement(i));
            }
        }
        else mp.put(this.INC_PLASTIC_MULTIPLIER, 0.0);
        return stress;
    }

    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#mountC(br.ufmg.dees.insane.analysisModel.AnalysisModel, br.ufmg.dees.insane.materialMedia.material.Material, java.util.HashMap)
     */
    public IMatrix mountC(AnalysisModel am, Material mat, HashMap hm) {
    	this.check(am,mat);
        return (am.getDualInternalVariableOperator(mat.getPs()));
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#check(br.ufmg.dees.insane.analysisModel.AnalysisModel)
     */
    public void check(AnalysisModel am) throws AnalysisModelNotValid {
        if (!(am instanceof PolarPlaneStrain ||
        	am instanceof PolarPlaneStress)	) {
            throw new AnalysisModelNotValid(this.getLabel()+" and "+am.getLabel());
        }
        
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#check(br.ufmg.dees.insane.analysisModel.AnalysisModel)
     */
    public void check(ContinuousPointModel am) throws ContinuousPointModelNotValid {
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#check(br.ufmg.dees.insane.materialMedia.material.Material)
     */
    public void check(Material mat) throws MaterialNotValid {
        if (!(mat instanceof PolarVonMises)) {
            throw new MaterialNotValid(this.getLabel()+" and "+mat.getClass());
        }
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#init(java.util.HashMap, java.util.HashMap, br.ufmg.dees.insane.analysisModel.AnalysisModel, br.ufmg.dees.insane.materialMedia.material.Material)
     */
    public void init(HashMap<Object, Object> previous, HashMap<Object, Object> cv, AnalysisModel am, Material mat) {
    	try {
			this.check(am);
			this.check(mat);
		} catch (AnalysisModelNotValid e) {
			e.printStackTrace();
		}catch (MaterialNotValid e) {
			e.printStackTrace();
		}
        previous.put(this.PLASTIC_DEFORM, new IVector(am.getNumberOfInternalVariables()));
        cv.put(this.PLASTIC_DEFORM, new IVector(am.getNumberOfInternalVariables()));
        cv.put(this.PLASTIC_MULTIPLIER, new Double(0.0));
        previous.put(this.PLASTIC_MULTIPLIER, new Double(0.0));
        cv.put(this.INC_PLASTIC_MULTIPLIER, new Double(0.0));
        cv.put(this.dFdSIG, new IVector(am.getNumberOfInternalVariables()));
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#init(java.util.HashMap, java.util.HashMap, br.ufmg.dees.insane.analysisModel.AnalysisModel, br.ufmg.dees.insane.materialMedia.material.Material)
     */
    public void init(HashMap previous, HashMap cv, ContinuousPointModel am, Material mat) {
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#update(java.util.HashMap, java.util.HashMap)
     */
    public void update(HashMap<Object, Object> previous, HashMap<Object, Object> cv) {
        previous.put(this.PLASTIC_DEFORM, cv.get(this.PLASTIC_DEFORM));
        previous.put(this.PLASTIC_MULTIPLIER, cv.get(this.PLASTIC_MULTIPLIER));
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#getNumberOfVariables()
     */
    public int getNumberOfVariables() {
        return 4;
    }
    
    /* (non-Javadoc)
     * @see br.ufmg.dees.insane.materialMedia.constitutiveModel.ConstitutiveModel#getNumberOfPreviousVariables()
     */
    public int getNumberOfPreviousVariables(){
        return 2;
    }
    
    /** Returns the yield function value from a given state of stress and internal variable. 
     *@return The yield function value.
     */
    private double yieldFunction(IVector s, AnalysisModel am, Material mat, HashMap mp){
     double J2 = 0.0;
     IMatrix P = am.getPMatrix();
       for(int i=0; i<am.getNumberOfInternalVariables(); i++)
           for(int j=0; j<am.getNumberOfInternalVariables(); j++)
               J2 += s.getElement(i) * P.getElement(i,j) * s.getElement(j);
       J2 *= 0.5;
       double f = Math.sqrt(3.0 * J2) - (((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).YIELDING).doubleValue() + ((PolarVonMises)mat).getMaterialValues(((PolarVonMises)mat).HARDENING_MODULUS).doubleValue()* ((Double)mp.get(this.PLASTIC_MULTIPLIER)).doubleValue()) ;
       return f;
    }
    
    /** Returns the gradient of yield surface w.r.t. the stress components. 
     *@return The gradient of yield surface w.r.t. the stress components.
     */
    private void yieldGrad(IVector s, AnalysisModel am, HashMap mp) { 
        double J2 = 0.0;
       IMatrix P = am.getPMatrix();
        for(int i=0; i<am.getNumberOfInternalVariables(); i++)
            for(int j=0; j<am.getNumberOfInternalVariables(); j++)
               J2 += s.getElement(i) * P.getElement(i,j) * s.getElement(j);
        J2 *= 0.5;
        double factor = 1.5 / Math.sqrt(3.0 * J2);
        for(int i=0; i<am.getNumberOfInternalVariables(); i++){
            ((IVector)mp.get(this.dFdSIG)).setElement(i, 0.0); 
            for(int j=0; j<am.getNumberOfInternalVariables(); j++)
                ((IVector)mp.get(this.dFdSIG)).setElement(i,((IVector)mp.get(this.dFdSIG)).getElement(i) + P.getElement(i,j) * s.getElement(j));
            ((IVector)mp.get(this.dFdSIG)).setElement(i,((IVector)mp.get(this.dFdSIG)).getElement(i) * factor);
        }
    }
    
    private void check(AnalysisModel am, Material mat){
    	try {
			this.check(am);
			this.check(mat);
		} catch (AnalysisModelNotValid e) {
			e.printStackTrace();
            System.exit(0);
		}catch (MaterialNotValid e) {
			e.printStackTrace();
            System.exit(0);
		}
    }
    
    public IMatrix mountC(ContinuousPointModel am, Material mat, HashMap cv) {
		// TODO Auto-generated method stub
		return null;
	}

	public IVector mountDualInternalVariableVector(IVector e, ContinuousPointModel am, Material mat, HashMap cv) {
		// TODO Auto-generated method stub
		return null;
	}

	public IMatrix mountCt(ContinuousPointModel am, Material mat, HashMap cv) {
		// TODO Auto-generated method stub
		return null;
	}

	public IMatrix mountCs(ContinuousPointModel am, Material mat, HashMap cv) {
		// TODO Auto-generated method stub
		return null;
	}

}
